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Abstract —Using fractional-calculus mathematics, a novel ap¬ 
proach is introduced to rewrite modified Bessel functions in 
series form using simple elementary functions. Then, a statistical 
characterization of the total receive-SNR at the destination, 
corresponding to the source-relay-destination and the source- 
destination link SNR, is provided for a general relaying scenario, 
in which the destination exploits a maximum ratio combining 
(MRC) receiver. Using the novel statistical model for the total 
receive SNR at the destination, accurate and simple analytical 
expressions for the outage probability, the bit error probability 
and the ergodic capacity are obtained. 


I. Introduction 


For an analytical investigation of AF cooperative commu¬ 
nication systems with a maximum-ratio-combining (MRC) at 
the destination, a statistical model of the total receive signal- 
to-noise ratio, SNR to t, is required, which is equivalent to the 
sum of SNRs corresponding to the Source-Destination (S-D) 
and the Source-Relay-Destination (S-R-D) link SNRs. To the 
best of our knowledge, a theoretical statistical model of SNR to t 
is yet unknown (but will be presented in this paper); numerous 
studies, however, consider the problem of finding a statistical 
model of SNR sr d- There are two main trends in the literature: 
either a single-branch system model is assumed, in which the 
receiver operates only on the relay transmission and, hence, no 
MRC is employed at the destination (e.g., see |TJ), or upper 
and/or lower performance bounds are considered assuming an 
MRC receiver at the destination (e.g., see ©)• 

To cope with the computations involving complicated math¬ 
ematical functions, a possible practical solution is to use an 
equivalent series representation of the functions (e.g. @). In 
fact, the appearance of the modified Bessel function of the 
second kind, K v {-), in the PDF of SNR sn j is the main source of 
trouble in AF-related calculations. We aim to substitute K u (-) 
with an equivalent series representation in this paper. However, 
as there is no appropriate and simple series representation 
of K v {-) available in the literature, motivated by Q, we 
derive a novel series representation of K v (-) in terms of 
simple elementary functions (such as x n e~ x ) using fractional 
calculus mathematics. With this result, the complex statistical 
model of SNR sr d turns out to be simple and easily tractable. 
Thereafter, the Cumulative Distribution Function (CDF) and 
the Probability Density Function (PDF) of SNR tot will be 
derived. 

The remainder of the paper is organized as follows: in 
Section [II] the system model is introduced. In Section III the 
fractional-calculus method is exploited to derive an equivalent 
series representation of K v (x). Based on the results derived 


in Section [III] the PDF and CDF of SNR tot are derived in 
Section |IV| and in Section [V] novel closed-form expressions 
are provided for some performance measures of variable-gain 
AF cooperative systems. Conclusion remarks are provided in 
Section [Vi] 


II. System Model 


We consider a two-hop variable-gain AF cooperative sys¬ 
tem. The source (S) sends data to the destination (D) by the 
help of an intermediate relay node (R). The destination “hears” 
both the source and the relay transmissions and employs a 
maximum-ratio-combining (MRC) receiver to jointly exploit 
all information available at the destination. Also, motivated 
by practical hardware constraints, it is assumed that the relay 
operates in half-duplex mode. 

The SNR at the output of the MRC receiver is 
SNR tot = 7 (|/i s d| 2 + |A src j| 2 ) = 'y\h eq \ 2 where 7 is the transmit 

\2 \h .12 _ IM 2 |M Q 


power at the source node. |/i sd | 2 and |A srd | 2 = 
represent equivalent S-D and S-R-D channel powers, respec¬ 
tively. |/i sd | 2 , | h sr | 2 and |/i rd | 2 are distributed exponentially 
with parameter A sd , A sr and A rd , respectively. The CDF of 
I A srd | 2 is 


^|^pW = l-2Ce- AsX ^ 1 (2C) 


(1) 


where ( = -\/Apx(x + I/7) and K v {-) the modified Bessel 
function of the second kind and v-\h order, Ap = A sr A rd and 
As = A sr + A rd ■ A proof of Q can be found in 0 and 0. 
The PDF of \h sld \ 2 is 


f\ hstd \2(x) = 2e~ XsX (\ P (2x + i)tf 0 (2C) + A s C^t(2C))- (2) 

The results in 0 and ([T]) do not, easily, lend themselves to 
further mathematical calculations (e.g. integration) as modified 
Bessel functions K v (■) appear. Hence, no statistical model for 
SNR tot is available in the literature so far. 

In the next section, an equivalent representation of K v (■) 
is derived that is based on a series-representation involving 
simple mathematical functions of the form x n e~ x . This novel 
equivalent representation of K v (■) is then used instead of the 
K{)(■) and Kj (■) functions appeared in 0 and 0 . It will be 
clear in the next sections that this approach paves the way for 
further theoretical analysis of AF relaying systems. 


III. Equivalent Series Representation of Modified 
Bessel Functions of Second Kind 
The mathematical concept of integration and differentiation 
of arbitrary (non-integer) order is called “fractional calculus”; 
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Fig. 1. Ki (( 3 x ) vs. series representation of Ki (( 3 x ) with k = 2 in FD 


foundations of the theory are discussed, e.g., in 0 . It will be 
used below to derive a simple novel equivalent representation 

of K v (fix). 

Theorem. Equivalent Representation of K v {fix) 

A modified Bessel function, K u (fix), of the second kind 
and v-th order, with v > 0 , can be represented by the infinite 
series 

oo n 

K v (fix) = ^2 ^2 A (^ *) ' {P x ) % ~ v i (3) 

n—0i —0 


with the coefficients 

(—1)* \J~f r(2z/) r(| + n — ^) L(n, i) 

A (is.n.l) = -:--j---^-j---;- 

2 u ~ l T(i - v) r(i + n + v) n\ 

that involve the Lah numbers (e.g., | 8 )) 

f n — 1\ n\ 
i — i ) a 


L{n , *) = 


for n,i> 0 , 


(4) 


(5) 


and the conventions L(0, 0) = 1; L(n, 0) = 0. 


Proof: Let s be a real non-negative number, i.e., s > 0 and 
s £ C. Let f{x) be continuous on x £ [0,oo) and integrable 
on any finite subinterval of x ^ 0. Then the Riemann-Liouville 
operator (e.g., |7j) of fractional integration is defined as 

I s {fix)} =-^ [ (x - t) s ~ l f{f)dt . (6) 

r («) Jo 

On the other hand, from [9] 3.471.4] we have 


[ {x-t) s - 1 t- 2s e~^ t dt = 
Jo 


T(s)/32 s -p 


e 2 = 


Ks ~^2x ) ' 


(7) 


Assuming f{t) = t~ 2s e~ the two integrals in (j 6 | and 
(|7]l are identical: this motivates the novel approach to derive 
an equivalent expression for K u (fix) by use of fractional 
integration. 

It follows from ([ 6 ]) and <[7]> that 


{s- 2 a e-^ x } = 


fi 2 s 


e^K s _i(f) 

7TX 2 AX 


( 8 ) 


The Leibniz rule for the Riemann-Liouville operator is 

I s {h{x)g{x)} = £ {h} D n {g} (9) 

n!I » 


where n is a non-negative integer, s + n is a non-negative 
fractional number and D n = By solving I s+n {h(x)} for 
h{x) = x~ 2s and D n {< 7 ( 2 ;)} for g{x) = e - ^/ 3 ^ the equivalent 
Bessel model |3]i will be derived. 

Let h{x) = x p , then according to (| 6 | 

I a {x p } = —*— / (a;-f ) a ~ 1 t p dt, (a > 0) (10) 

r(a) Jo 


r(i+p) 


r p+a 


rn , ^ x- • ( 11 ) 

r(l +p + a) 

with r(-) the Gamma-function. Suppose that p = — 2s and 
a = s + n, then 


i(s+u) | x ~ 2s | = JXi—^ 
1 J T(1 — s + n) 


( 12 ) 


With g(x) = e P/ x in (J9J), D n {e P/ x } can be found to 
equal 


D 7 


{e-4/*} 


/I \n n 

= e -/»/*L_2_£(_i y L {n,i){/3/ x y, (13) 


X" 


i =0 


with L(n,i ) defined in ([5]>; this result is taken from 1101. 

By substituting ( [T2] i and ( p~3] > into <[ 8 ]i and ([9]) it is straight¬ 
forward to obtain 




e -4 K ,A = 

•—X2x' Z- n! r(s)r(l + n - s) 

n 


Xe -P/ X LJL J2i-iyL{n, i)(fi/x)\ (14) 


X" 


i =0 


Changing the variable x —> 2, assuming 1 — 2s = 2v and 
exploiting A'_„ = K v , the result is the infinite series in (|3j. 


It should be made clear that the above representation of 
K„(fix) is not valid for v = |0, 4, |, • • • }. That is because 


T(2v) and r(| + ro — v) in (|4j diverge to ±oo. However, 
for the case of v = 0, one can compute K 0 (fix ) using 
the equivalent representation of K-\ [fix) and K 2 (fix ) by 


K v (fix) = K u _ 2 (fix) + K v -i{Px) 

from (m 10.38.4], 


that is obtained 


Finite Series Representation of K u (fix) 

The equivalent representation of K u (fix) may significantly 
simplify computations involving K u (fix), as the series in ([3]) 
contains the variable x only in the simple function-template 

x i-v e -Px 

that can, e.g., be easily integrated. The series 
representation contains, however, an infinite number of terms 
that can’t be computed in practical applications. 

Fortunately, the series representation of K„ (fix) is rather 
accurate for a finite number of terms as defined as follows: 

k n 

K v (fix) = e-? x ^2 E A ^’ n ’ 0 • + e ( 15 > 

n—0 2—0 
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with 


e = e-^ ^ ^ A(^, n, i) • (/3x) i_ ". (16) 

n=fc+l i=0 

The first term on the right-hand side of ( p~5| ) represents the 
actual function to approximate K u (/3x), and e represents 
the truncation error that can be neglected. Fig. [T] illustrates 
numerical values of the finite series representation of K\ (ftx) 
(with k = 2 in ( fl5] l) for various values of ft (dashed lines) 
and also the theoretical fully accurate values of Ki(ftx) (solid 
lines). It is clear from the figure that the finite series for 
K„(/3x) with only k = 2 merges with theoretical K v (ftx) 
with high accuracy. 

Convenient Representation for Practical Use of the Trun¬ 
cated Series: For practical use, it is convenient to re-write 


(15 i such that one of the sum-operators is included in the 
series coefficients. For this, the inner sum over i in © is 
evaluated row-wise, with the sum-index n counting the rows. 
This structure is then summed up column-wise and the result 
can be written as 


Kft(fix) = 




F E (E A ^4,<7)) •(/? 


e (17) 


q —0 l—q 


— ^iy,k,q 

As long as k is limited, the above series representation ( fT7| ) 
can always be used to replace 0 without any convergence 
problems but still the truncation error e can be made arbitrarily 
small. Numerical values for the coefficients a v ^,q are given in 
Table [I] for the first-order (v = 1) modified Bessel function of 
the second kind ITi(-); it should be noted that the coefficient 
with g-index 1 is always found to equal ai^,i = 2 k+i ■ 


IV. Distribution of Equivalent Channel Power 

Assuming an MRC receiver at the destination, the total 
receive SNR is the sum of SNRs corresponding to the S- 
D and the S-R-D links, i.e., SNR tot = SNR st j + SNR sli j = 
7 (|^sd| 2 + |^srd| 2 )- The equivalent channel power-gain can be 
written as \h eq \ 2 = |/i sr d| 2 + |/i s d| 2 . and so, for the CDF of 
|F eq | 2 we obtain 

F\h eq \ 2 ( x ) = P (IM 2 + I*srd| 2 < X ) 


— A s de 


— A sd a: 


* XsdU -F\ h M u )du, (18) 


with F|, isrd | 2 (.) derived in 0 For simplicity we will restrict 
calculations below to the high “transmit-SNR” regime but it 
will be demonstrated, by simulations (e.g., Fig. [3]) that this is 
justified because it leads to very accurate numerical results, 
even in the low-SNR region. 

By plugging ([l} into (18 1 , and assuming “high SNR”, (181 
simplifies to 


iKl*(r) = 1 - e- A " r - 2A sd yA^e- Asdr 


(19) 


[ xe- {Xst+Xti ~ Xsd)x K 1 (2^K\dx)dx. 

Jo 


The integral in ( [T9] ) is non-trivial and does not seem to have 



Fig. 2. PDF f\heq\ 2 ( x ) us i n g the theoretical model derived in |24|; Monte 
Carlo simulations and the bounding technique proposed in theiiterature 
(e.g. jT2j), all with A sr = A rd = 1. 


closed-form solution. However, using the results from Section 
m the integral can be rewritten as follows 

V = [ 2yJ A sr A rd ue _(Ast+Ard_Asd) "/\i( 2yJX SI X ld u)dx 

Jo 

(2-s/A S rArd) 9 (?! Cll,k,q 
(Asrd - A sd )9 +1 

x (i - jr A ri -A sd )x) (20) 

c=0 C ‘ 

where A sr d = (v%^ + -\/Ajd) 2 - The second step is obtained 
by using the series representation of Ki(2\fXfX U ix) derived 
in As the series is truncated (for k limited), this is an 
approximation, indicated by the use of instead of strict 
equality; the truncation error can, however, be made arbitrarily 
small by choosing proper k. By substituting ( [20] ) into m k 
is straightforward to obtain Fj^p (x) as 

k q 

F\ heq \2(x) « 1 - Ae~ XsiX + Y.Y, ^ C e _Asrd " ( 21 ) 

q—0 c—0 

where coefficients A and B are independent of x, defined as 


A 

B 


g=0 


Asd(2"\/A sr A r d) 9 g! 

(Asrd — A sd )' ?+1 


Asd(2"\/ A sr A r d) 9 <z! ti\ k,q 

C! (A srd - A sd )9- C+1 


( 22 ) 

(23) 


The PDF of |F e q| 2 is the derivative of Fj^p (x) in (21 1 w.r.t x, 
which is easy to calculate as the series representation involves 
simple elementary functions only: 


k q 

/|N q | 2 (•*-) ~ -4A s de- Asda: + EE B(cx c 1 — A srd a; c )e Asrda: . 

q—0 c—0 

(24) 

Fig. |2jillustrates the accuracy of f\h cq \ 2 (x) using the expression 
derived in (|24j) (solid line) by a comparison with histogram- 
results obtained from Monte Carlo simulations. 
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a l,fc,0 

ai.fc.l 

0-1,k,2 

a l,fc,3 

0-1,k,A 

a l,k,5 

a l,k,6 

«l,fc,7 

O-l,k,8 

a l,k,9 

Ol,k,10 

k = 2 

1 

4/5 

-0.1333 



II 

Cn 

1 

10/11 

-0.4237 

0.1824 

-0.0375 

2.693 

xlO -3 


O 

II 

1 

20/21 

-0.7047 

0.7239 

-0.5000 

0.2111 

-5.415 

xlO " 2 

8.375 

xlO " 3 

-7.55 

xlO " 4 

3.619 

xlO " 5 

-7.0724 

xlO -7 


Coefficients a Ut k, q in O 


FOR V = 1 WITH FOUR DIGITS OF ACCURACY 



Fig. 3. Ergodic Capacity (nats/sec/Hz) of the system with multiple relays 
and equal transmit power from the relays. 


In the literature (e.g. 1121), when considering relaying 


systems with MRC receiver at the destination, a method is 
proposed for estimating the statistics of SNR SI( i based on 
the bounding technique SNR srd = mm(SNR sr , SNR rl ]). Conse¬ 
quently, SNR sr( i has as an exponential distribution with param¬ 
eter Asr+Ard- Note that, although this method greatly simplifies 
the calculations by avoiding modified Bessel functions in the 
formulations, accuracy is sacrificed: Fig. [2] (dashed line) shows 
the result for f\ heq ^(x) when using the bounding method as 
well. It is clear that the approach presented in this paper results 
in more accurate solutions in comparison with traditional 
bounding techniques. 


V. Performance Analysis 
The CDF in 0 corresponds to the outage probability. 
Due to lack of space no illustration is provided for outage 
probability but with the accuracy of the PDF in Fig. [2] the 
accuracy of outage probability is ensured as well. 

bit error probability (BEP) 


24 


Given the PDF of |h . eq | 2 in 
can be obtained as p b = | f Q eiici^Nyx) f\h eq \ 2 (x)dx■ From 


113 7.1.19], f eric(^/jx) = — 1//2 e TX . Then, using 


integration by parts, pb can be written according to 




7 - 4 2 


7 + A s 


sd 


1 — k q 

v q= 0 c—0 

BT (c+i) 

)■ (25) 

(7 + Asrd) C+ 5 

of an AF 

cooperative 

system 


can be obtained using the definition of capacity as 


C av = 5/o°° lo g( 1 +l x )f\h^\ 2 (x)dx, that can easily be 
solved in closed-form by substituting /|^ e | 2 (x) from m 
and exploiting [9] 4.222.8]. Due to space limit, the closed 
form expression is not provided in this letter; however, Figj3] 
illustrates the Ergodic capacity obtained from the proposed 
approach. Also, the proposed approach has been examined 
for multi relay scenarios as well. The excellent agreement 
between Monte-Carlo Simulations and closed-form solution 
proves the correctness of the results. 

VI. Conclusions 

Based on fractional calculus mathematics, a novel power 
series representation of K v (■) is introduced that has simple 
elementary functions of the form x n e~ x as its basic com¬ 
ponents. Based on that the performance of AF cooperative 
systems is investigated and accurate results are obtained. 
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